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ABSTRACT 

This  paper  describes  novel  experimental  methods  aimed  at  understanding  the  fundamental  phenomena 
governing  the  motion  of  lightweight  vehicles  on  dry,  granular  soils.  A  single-wheel  test  rig  is  used  to  empirically 
investigate  wheel  motion  under  controlled  wheel  slip  and  loading  conditions  on  sandy,  dry  soil  Test  conditions 
can  be  designed  to  replicate  typical  field  scenarios  for  lightweight  robots,  while  key  operational  parameters 
such  as  drawbar  force,  torque,  and  sinkage  are  measured.  This  test  rig  enables  imposition  of  velocities,  or 
application  of  loads,  to  interchangeable  running  gears  within  a  confined  soil  bin  of  dimensions  1.5  m  long,  0.7  m 
wide,  and  0.4  m  deep.  This  allows  testing  of  small-scale  wheels,  tracks,  and  cone  or  plate  penetrators.  Aside 
from  standard  wheel  experiments  (i.e.,  measurements  of  drawbar  force,  applied  torque,  and  sinkage  during 
controlled  slip  runs)  two  additional  experimental  methodologies  have  been  developed.  The  first  relies  on  high¬ 
speed  imaging  of  the  wheel-soil  interface  and  the  use  of  particle  image  velocimetry  (PIV)  to  measure  micro¬ 
scale  terrain  kinematics.  The  second  experimental  methodology  consists  of  a  custom  force  sensor  array  located 
at  the  wheel-terrain  interface.  The  sensors  allow  explicit  measurement  of  normal  and  shear  forces  (and, 
therefore,  estimation  of  normal  and  shear  stresses)  at  numerous  discrete  points  along  the  wheel-soil  interface. 
Experimental  measurements  gathered  by  these  test  methodologies  are  to  be  compared  against  well-established 
semi-empirical  models,  to  validate  and  understand  limitations  of  the  models  and  propose  improvements. 


INTRODUCTION 

Robotic  vehicles  are  frequently  deployed  in  unwelcoming, 
hazardous  environments.  From  military  robots  to  planetary 
rovers,  vehicle  mobility  is  a  key  aspect  of  mission  success. 
Several  models  for  traction  modeling  of  tracked  and  wheeled 
vehicles  have  been  developed  in  the  past  decades;  however, 
a  comprehensive  understanding  of  soil  behavior  under 
running  gear  is  still  missing  to  date.  The  work  of  Bekker  and 
Wong,  which  began  in  the  1950’s,  has  laid  the  foundation 
for  modern  terramechanics.  The  application  of  classical 
results  from  plasticity  theory,  combined  with  semi-empirical 
formulations,  has  provided  satisfactory  solutions  to  the 
problem  of  mobility  modeling  for  large,  heavy  vehicles. 
However,  the  expanded  use  of  lightweight  vehicles 


(especially  man-portable  robotic  vehicles)  has  called  for  a 
new  effort  in  modeling  vehicle-terrain  interaction  problems. 
In  fact,  some  researchers  have  suggested  that  classical 
models  are  of  questionable  utility  when  applied  to  vehicles 
one  order  (or  more)  of  magnitude  smaller  than  tanks, 
Humvees,  large  trucks,  and  the  like  [1]. 

This  paper  will  describe  novel  experimental  methods 
aimed  at  understanding  the  fundamental  phenomena 
governing  the  motion  of  lightweight  vehicles  on  dry, 
granular  soils.  A  single-wheel  test  rig  is  used  to  empirically 
investigate  wheel  motion  under  controlled  wheel  slip  and 
loading  conditions  on  a  sandy,  dry  soil  (Figure  1).  Test 
conditions  can  be  designed  to  replicate  typical  field 
scenarios  for  lightweight  robots,  while  key  operational 
parameters  such  as  drawbar  force,  torque,  and  sinkage  are 
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measured.  This  test  rig  enables  imposition  of  velocities,  or 
application  of  loads,  to  interchangeable  running  gears  within 
a  confined  soil  bin  of  dimensions  1.5  m  long,  0.7  m  wide, 
and  0.4  m  deep.  This  allows  testing  of  small-scale  wheels, 
tracks,  and  cone  or  plate  penetrators. 

The  soil  under  investigation  has  been  fully  characterized 
with  a  series  of  direct  shear  tests  (ASTM  D3080)  and 
penetration  tests.  Direct  shear  tests  were  performed  to 
estimate  soil  shearing  parameters  such  as  cohesion,  angle  of 
internal  friction,  and  shear  modulus.  Penetration  tests, 
although  not  standard  tests,  were  performed  to  evaluate 
‘Bekker’  parameters,  necessary  for  characterization  of 
pressure-sinkage  behavior  of  the  soil  under  the  methodology 
described  by  Wong  [2]. 

The  aforementioned  experiments  represent  a  typical 
experimental  approach  to  macro-scale  characterization  of 
wheel-soil  interaction.  However,  the  application  of  classical 
terramechanics  model  to  lightweight  vehicles  may 
potentially  show  discrepancy  between  experiments  and 
predictions,  warranting  the  development  of  new  methods  to 
probe  the  fundamental  mechanics  of  a  small  robot’s 
interaction  with  soil. 

To  this  end,  two  additional  experimental  methodologies 
have  been  developed.  The  first  relies  on  high-speed  imaging 
of  the  wheel- soil  interface  and  the  use  of  particle  image 
velocimetry  (PIV)  to  measure  micro-scale  terrain 
displacement  (Figure  1).  This  methodology,  although 
confined  to  plane  strain  cases,  allows  measurement  of  soil 
flow  velocities,  and  observation  of  the  formation  of  shear 
bands  beneath  the  wheel/track.  Though  this  method  does  not 
explicitly  permit  calculation  of  the  velocities  of  individual 
soil  particles,  it  does  allow  estimation  of  a  regularly- spaced 
velocity  field  in  the  soil.  While  such  visualization 
techniques  have  been  widely  employed  in  the  field  of 
experimental  fluid  mechanics,  their  application  to  the  study 
of  soils  is  a  relatively  new  development  [3,  4]. 

The  second  experimental  methodology  is  intended  to 
complement  the  PIV-based  soil  kinematics  analysis.  It 
employs  a  custom  force  sensor  array  located  at  the  wheel- 
terrain  interface.  The  force  sensors  are  strain  gauge-based 
flexural  elements  with  interchangeable  interface  surfaces 
that  are  designed  for  integration  with  wheels  or  other 
running  gear.  The  sensors  allow  explicit  measurement  of 
normal  and  shear  forces  (and,  therefore,  estimation  of 
normal  and  shear  stresses)  at  numerous  discrete  points  along 
the  wheel-soil  interface.  When  coupled  with  PIV-derived 
kinematic  data,  this  allows  for  a  richer  characterization  of 
soil  loading  and  failure  regimes  than  would  be  possible  with 
either  kinematic  or  pressure  information  alone.  In  particular, 
this  experimental  methodology  allows  joint  visualization  of 
the  soil  displacement  in  the  bulk  soil  medium,  and 
measurement  of  shear  and  normal  stress  at  points  along  the 
interface.  This  could  lead  to  development  and  validation  of 


novel  constitutive  relations  describing  soil  behavior  under 
loading  imposed  by  running  gear. 

Experiments  have  shown  that  soil  failure,  at  certain  slip 
levels,  is  qualitatively  different  under  cases  of  low  vertical 
load  (which  is  typical  for  lightweight  robots)  compared  to 
cases  of  high  vertical  load  (typical  for  large  ground 
vehicles).  Also,  soil  flow  patterns  have  been  observed  to 
exhibit  periodic  failure  phenomena,  giving  rise  to  interesting 
features  such  as  surface  ripple  formation.  These  results, 
obtained  through  PIV  analysis,  provide  deeper 
understanding  of  the  mechanics  of  traction  generation. 
Experimental  measurements  gathered  by  these  test 
methodologies  are  compared  against  the  results  from  well- 
established  semi-empirical  models,  to  understand  limitations 
of  these  models  and  propose  modifications  and 
improvements. 


Figure  1:  CAD  drawing  of  the  terramechanics  testbed  showing 
the  imager  for  PIV  experiments  (top).  Actual  PIV  setup  with  the 
high  speed  camera  and  two  flood  lights  (bottom). 
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SINGLE  WHEEL  TESTBED  DESCRIPTION 

The  Robotic  Mobility  Group  at  MIT  has  designed  and 
fabricated  a  multipurpose  terramechanics  rig  based  on  the 
standard  design  described  by  Iagnemma  [5].  The  testbed  is 
pictured  in  Figure  1  and  is  composed  of  a  Lexan  soil  bin 
surrounded  by  an  aluminum  frame  where  all  the  moving 
parts,  actuators  and  sensors  are  attached.  A  carriage  slides  on 
two  low- friction  rails  to  allow  longitudinal  translation  while 
the  wheel  or  track,  attached  to  the  carriage,  is  able  to  rotate 
at  a  desired  angular  velocity.  The  wheel  mount  is  also  able 
to  translate  in  the  vertical  direction.  This  typical  setup  allows 
control  of  slip  and  vertical  load  by  modifying  the 
translational  velocity  of  the  carriage,  angular  velocity  of  the 
wheel,  and  applied  load.  Horizontal  carriage  displacement  is 
controlled  through  a  toothed  belt,  actuated  by  a  90W  Maxon 
DC  motor  while  the  wheel  is  directly  driven  by  another 
Maxon  DC  motor.  The  motors  are  controlled  thorough  two 
identical  Maxon  ADS  50/10  4-Q-Dc  servoamplifiers.  The 
carriage  horizontal  displacement  is  monitored  with  a  Micro 
Epsilon  WPS-1250-MK46  draw  wire  encoder  while  wheel 
vertical  displacement  (i.e.,  sinkage)  is  measured  with  a 
Turck  A50  draw  wire  encoder.  A  6-axis  force  torque  ATI 
Omega  85  transducer  is  mounted  between  the  wheel  mount 
and  the  carriage  in  order  to  measure  vertical  load  and 
traction  generated  by  the  wheel.  Finally,  a  flange -to-flange 
reaction  torque  sensor  from  Futek  (TFF500)  is  used  to 
measure  driving  torque  applied  to  the  wheel.  Control  and 
measurement  signals  are  handled  by  a  NI  PCIe-6363  card 
through  Lab  view  software. 

The  rig  is  capable  of  approximately  1  meter  of  horizontal 
displacement  at  a  maximum  velocity  of  approximately  120 
mm/s  with  a  maximal  wheel  angular  velocity  of 
approximately  40  deg/s.  The  bin  width  is  0.6  meters  while 
the  soil  depth  is  0.16  meters.  Considering  the  wheel  sizes 
and  vertical  loads  under  study,  these  physical  dimensions  are 
sufficient  for  eliminating  boundary  effects.  Moreover,  the 
same  testbed,  with  some  adaptations,  can  be  used  to  perform 
soil  penetration  tests  and  analyze  different  running  gears 
(e.g.,  both  wheels  and  tracks). 

For  the  experiments  described  in  this  paper,  the  Mojave 
Martian  Simulant  (MMS)  was  employed  as  a  test  medium 
[6].  MMS  is  a  mixture  of  finely  crushed  and  sorted  granular 
basalt  intended  to  mimic,  both  at  chemical  and  mechanical 
level,  Mars  soil  characteristics.  MMS  particle  size 
distribution  spans  from  micron  level  to  mm  level  with  80% 
of  particles  above  the  1 0  micron  threshold. 


GRANULAR  SOIL  PARTICLE  IMAGE 
VELOCIMETRY 

Particle  image  velocimetry  (PIV)  describes  an 
experimental  method,  based  on  image  cross-correlation 
techniques,  used  for  the  determination  of  flow  velocity 
fields.  The  use  of  PIV  for  the  calculation  of  fluid  velocities 
initially  emerged  in  the  19804s  [7,  8].  Since  then,  PIV  has 
played  an  important  role  in  many  fluid  mechanics 
investigations  [9].  Two  of  the  main  advantages  of  PIV  over 
other  methods  for  the  measurement  of  velocity  (e.g.  hot- 
wire-velocimetry,  Pitot  tubes  etc.)  are  that  it  is  non-intrusive, 
and  allows  for  relatively  high  resolution  measurements  over 
an  extended  spatial  domain. 

During  fluid-based  PIV  analysis,  the  fluid  is  typically 
seeded  with  marker  particles  that  refract,  absorb,  or  scatter 
light,  have  a  high  contrast  with  the  fluid,  and  do  not  interrupt 
the  fluid  flow.  Imaging  is  performed  at  high  speed  over  an 
area  of  the  flow  illuminated  by  a  light  source,  typically  a 
pulsed  laser.  Captured  images  are  post-processed  with 
algorithms  that  perform  frame-to-frame  feature  tracking  and 
calculation  of  flow  velocity  fields. 

PIV  is  also  a  useful  method  for  measuring  soil  motion, 
with  the  notable  constraint  that  soil  is  typically  observed 
through  a  glass  sheet,  limiting  the  resulting  analysis  to  plane 
strain  scenarios.  The  natural  granular  texture  of  soils  often 
generates  an  intensity  pattern  that  can  be  readily  traced  by 
PIV-algorithms,  without  the  use  of  marker  particles.  Also, 
incandescent  light  can  generally  be  used  for  illumination. 

Granular  PIV  has  recently  been  employed  in  several 
applications,  including  the  analysis  of  grains  in  converging 
hoppers  [10],  study  of  flowing  granular  layers  in  rotating 
tumblers  [11],  investigation  of  granular  avalanches  [12], 
analysis  of  soil  motion  caused  by  the  movement  of  animals 
[13],  the  study  of  burrowing  behavior  of  razor  clams  [3],  and 
in  the  study  of  wheel-soil  interaction  [4,  14].  The  analysis  of 
soil  motion  beneath  a  driven  wheel  via  quantitative  analysis 
of  successive  temporal  images  was  first  introduced  by  Wong 
[15].  However,  the  experimental  capabilities  of  that  study 
did  not  allow  for  high-speed  image  capture,  limiting  the 
accuracy  and  practical  utility  of  the  method. 

Soil  motion  analysis  can  be  broken  down  into  four  main 
steps:  1)  image  acquisition,  2)  image  pre-processing,  3) 
image  cross-correlation  (PIV),  and  4)  velocity  field  post¬ 
processing.  These  steps  are  briefly  described  here,  and 
methods  for  parameter  selection  are  presented.  Note  that  in 
the  following,  the  Matlab-based  PIVlab  software  is 
employed  [16]. 
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Figure  3:  Two  examples  of  image  canonical  transformations 
used  to  evaluate  PIV  settings.  Nine  image  transformations  for 
coarse  and  fine  soil  textures  were  used  to  evaluate  PIV 
accuracy. 


PIV  IMAGER  CONFIGURATION 

The  accuracy  of  PIV  strongly  depends  on  the  quality  of  the 
captured  images.  For  these  experiments  the  testbed  was 
fitted  with  a  2.54  cm  thick  tempered  glass  wall  while  the 
running  gear  was  operated  flush  against  this  surface  (see 
Figure  1).  Both  wheels  and  tracks  have  been  analyzed  with 
this  testbed,  however  this  paper  describes  results  from  rigid 
wheel  testing. 

Image  sets  for  the  PIV  measurement  were  captured  with  a 
Phantom  7  high-speed  camera.  The  Phantom  7  is  able  to 
record  grayscale  images  at  the  maximum  resolution  of 
800x600  pixels  at  a  maximum  frame  rate  of  6688  fps.  The 
camera  was  placed  perpendicular  to  the  front  glass  wall  (see 
Figure  1)  at  a  distance  of  52  cm,  while  the  focal  length  was 
set  to  77  mm  (a  zoom  lens  was  used)  resulting  in  an  image 
capture  region  of  approximately  15  x  1 1.25  cm.  It  should  be 
noted  that  determination  of  image  capture  region  size  is 
largely  dictated  by  the  particular  experimental  conditions. 
Here,  the  image  capture  region  was  chosen  in  order  to 
conservatively  bound  the  region  of  soil  that  would  undergo 
motion  when  subjected  to  wheel  passage  on  the  soil  surface. 
Two  250W  Lowel  Pro-Light  photography  flood  lights  were 
placed  on  either  side  of  the  camera  at  an  angle  of  45° 
towards  the  object  plane,  and  provided  approximately 
homogeneous  illumination  of  the  soil.  By  using  two  laterally 
positioned  light  sources,  reflections  and  shadows  can  be 
significantly  diminished. 


PIV  IMAGE  PREPROCESSING 

The  performance  of  PIV  cross-correlation  algorithms 
generally  improves  when  images  are  of  high  contrast,  feature 
dense,  and  have  low  noise.  In  practice,  images  are  subject  to 
nonuniform  illumination,  image  sensor  noise,  and  lack  of 
natural  contrast  in  the  granular  material,  all  of  which  can 
degrade  PIV  algorithm  performance.  Various  image  pre¬ 
processing  methods  were  investigated  to  understand  their 
effect  on  algorithm  performance.  These  include  commonly- 
employed  algorithms  such  as  contrast  limited  adaptive 
histogram  equalization,  high  pass  filtering,  and  clipping  and 
intensity  capping. 

To  systematically  investigate  the  effect  of  these 
preprocessing  methods  on  PIV  algorithm  performance,  test 
image  segments  of  the  Mars  regolith  simulant  with 
dimensions  256  x  256  pixels  were  captured,  then 
synthetically  deformed  in  canonical  directions.  Since  the 
particle  distribution  in  the  soil  under  investigation  is  locally 
inhomogeneous,  two  distinct  image  segments  were  captured 
in  order  to  adequately  represent  typical  apparent  grain 
distributions  in  the  MMS  simulant.  This  resulted  in  one 
image  populated  by  relatively  large  grains  and  one  populated 
by  relatively  small  grains  (Figure  2).  Synthetic  deformation 
of  the  image  was  performed  as  a  means  of  generating  a 
ground  truth  for  cases  of  linear  translation  (1-4  pixels  in  both 
horizontal,  vertical,  and  diagonal  directions),  rotation  (1-8 
degrees  in  clockwise  and  counter-clockwise  directions), 
shear  (1-4  pixels  of  relative  motion  between  upper  and  lower 
image  halves),  and  simple  shear  (1-4  pixels  of  motion  of 
upper  edge  of  image)  (Figure  3).  Since  the  pixel  shift  for 
each  deformation  was  controlled,  this  methodology  allowed 
quantitative  evaluation  of  PIV  algorithm  results.  An  error 
metric  was  computed  by  computing  the  average  difference, 
over  all  points  in  the  PIV  velocity  field,  between  the  velocity 
vector  calculated  through  PIV  and  the  true  velocity  vectors. 

PIV  IMAGE  CROSS-CORRELATION 

In  PIV,  images  are  divided  into  small  interrogation 
windows  (IW)  and  then  analyzed  to  compute  the  probable 
displacement  between  successive  images  for  each  IW  using 
cross-correlation  techniques.  This  results  in  an  equally 
spaced  field  of  calculated  velocity  vectors.  The  probable 
displacement  is  determined  by  using  the  cross-correlation 
function: 

Ri,'(x,y)  =  £;=_/, /(i,;)  /'(*  +  x,j  +  y )  (1) 

where  /  is  the  intensity  of  the  first  image  and  /’  the  intensity 
of  the  second  image.  A  detailed  description  of  PIV  theory 
can  be  found  in  [17].  Particle  density,  image  resolution,  and 
IW  size  are  interconnected  parameters  that  must  be  carefully 
selected  to  optimize  performance.  Based  on  experimental 
investigations,  Keane  and  Adrian  [18]  defined  empirical 
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Figure  4:  Comparison  of  velocity  calculated  through  PIV  and 
measured  with  a  draw  wire  encoder. 


Figure  5:  Soil  trajectories  calculated  from  velocity  field  obtained 
through  PIV  analysis.  Visual  inspection  showed  that  PIV  yielded 
tracking  of  soil  regions  on  the  order  of  0.5-1  mm  after  translations  of 
several  centimeters. 


rules  for  optimal  PIV  setup.  The  reader  is  referred  to  the 
above  paper  for  more  details.  For  the  results  presented  here, 
the  following  settings  were  employed:  25  fps,  final  IW  size 
of  16,  CLAHE  filtering  with  kernel  size  of  40  pixels.  A 
more  complete  description  of  the  PIV  settings  and  analysis  is 
presented  in  [19]. 

VELOCITY  FIELD  POST-PROCESSING 

The  raw  velocity  field  produced  by  PIV  calculations  can 
contain  spurious  vectors  (outliers).  These  outliers  can  be 
caused  by  noise,  inappropriate  interrogation  settings,  and 
accidentally  matched  patterns.  Hence,  to  improve  results, 
rejection  of  these  outliers  and  interpolation  of  missing  data 
points  can  be  performed  in  a  post-processing  stage  through 


filtering.  Filters  for  the  rejection  of  outliers  can  primarily  be 
divided  into  two  separate  classes:  global  and  local  methods. 
Global  filters  commonly  employ  a  simple  thresholding 
method,  with  the  threshold  value  selected  by  an  operator 
with  empirical  or  theoretical  domain  knowledge.  If  elements 
of  the  velocity  field  exceed  the  threshold,  this  element  is 
removed  from  the  results.  Local  filters  are  primarily  based 
on  relative  differences  between  surrounding  vectors,  rather 
than  absolute  values.  A  local  filter  calculates  the  mean  and 
standard  deviation  of  the  velocity  for  a  selected  kernel  size 
around  each  vector.  If  the  velocity  exceeds  certain 
thresholds,  the  vector  is  rejected.  For  the  results  presented 
here,  a  5x5  kernel  with  a  threshold  of  8  times  the  standard 
deviation  was  used  for  post-processing. 

VALIDATION  AND  VERIFICATION 

The  synthetically  deformed  image  was  determined  to  be  a 
useful  ground  truth  for  determining  appropriate  PIV 
operational  parameters.  However,  validation  of  the  PIV 
algorithm  performance  was  also  pursued  on  two  sets  of  test 
data  that  were  physically  relevant  to  the  running  gear-soil 
interaction  case. 

The  first  test  consisted  of  calculating  the  velocity  via  PIV 
of  a  2.5  cm  thick  steel  plate  performing  a  soil  penetration 
test.  The  ground  truth  velocity  of  the  plate  was  externally 
measured  by  numerically  differentiating  the  output  of  the 
draw  wire  encoder  (which  nominally  provides  a  position 
measurement).  To  obtain  a  plate  velocity  measure  from  PIV, 
an  average  of  the  velocities  was  computed  over  a  rectangular 
region  of  interest  aligned  with  the  moving  plate. 

Figure  4  shows  a  comparison  of  the  plate  velocity  as 
determined  from  PIV  calculations  and  the  velocity  measured 
by  the  draw  wire  encoder.  The  average  percent  error  (for  the 
best  settings)  between  these  measurements  was  below  1%.  It 
should  be  noted  that,  for  this  test  case,  the  PIV  algorithm  is 
not  performing  calculations  on  the  granular  soil,  but  rather 
the  steel  plate  edge.  However,  this  test  remains  of  interest 
since  the  soil  in  contact  with  the  plate  necessarily  moves  at 
the  same  velocity. 

The  second  test  consisted  of  calculating  the  time  evolution 
of  motions  of  discrete  features  associated  with  MMS 
simulant  soil  beneath  a  driven  rigid  wheel.  Trajectories  s(t) 
are  calculated  for  a  grid  of  9  x  6  regions  of  interest  over  the 
soil  area.  The  time  evolution  of  the  positions  of  the  center  of 
the  regions  of  interest  was  computed  by  integrating  the 
velocities  with  a  fourth  order  Runge-Kutta  method. 


s(t)  =  I  v(t)dt  (2) 

Jo 
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Figure  6:  Working  scheme  of  the  custom  force  sensor  for  interfacial  stress  measurement  (top  left).  Five  sensors  are  distributed  from  the  wheel 
median  axis  to  the  wheel  edge  (bottom  left).  Sensors  are  rigidly  connected  to  the  wheel  hub  (right) 


The  motion  of  these  tracked  regions  were  compared  to 
trajectories  of  individual  soil  particles  that  are  large  enough 
to  be  manually  tracked  from  frame  to  frame,  thereby 
providing  a  qualitative  performance  evaluation.  Also,  the 
calculation  of  feature  trajectories  is  useful  for  illustrating 
soil  flow  when  subjected  to  various  loading  conditions. 

Figure  5  displays  the  trajectories  computed  while  the 
wheel  was  advancing  at  17  deg/s  with  30%  slip.  Note  that 
the  area  above  the  soil  surface  was  masked  during  pre¬ 
processing,  and  hence  these  features  remain  at  their  original 
location.  The  squares  show  the  final  position  of  the  tracked 
features  and  the  lines  represent  the  motion  evolution. 
Manual  inspection  showed  that  the  selected  PIV  yielded 
tracking  of  soil  regions  in  the  order  of  1-2  pixels, 
corresponding  to  0.5-1  mm  after  translations  of  several 
centimeters. 

WHEEL-TERRAIN  INTERFACE  FORCE  SENSOR 
DESCRIPTION 

Measurement  of  the  normal  and  shear  stress  acting  on  a 
moving  wheel  is  important  for  empirical  testing  and 
validation  of  models  describing  interfacial  phenomena. 
While  numerous  COTS  sensors  exist  for  measuring  pressure 
[20],  the  authors  are  unaware  of  any  available  sensors  that 
can  measure  both  pressure  and  shear  stress,  at  a  scale  and 
resolution  suitable  for  investigation  of  the  interaction 
mechanics  of  small,  lightweight  vehicle  running  gear  and 
deformable  soil. 

Therefore,  a  custom  sensor  array  was  designed  and 
fabricated  (Figure  6).  Each  sensor  is  a  solid-state  L-shaped 
aluminum  flexure  instrumented  with  two  full  bridge  strain 


gages.  The  sensor  is  mounted  rigidly  to  the  running  gear, 
and  its  interface  element  is  exposed  to  the  soil.  The 
interface  element  is  generally  subjected  to  normal  (TV)  and 
shear  ( T)  loading.  These  forces  cause  the  flexure  elements 
to  deflect  in  a  linear  elastic  manner.  From  measured 
deflection,  and  given  prior  calibration  data,  the  applied 
forces  can  be  uniquely  computed.  (Axial  strain  is 
intrinsically  rejected  by  the  full  bridge  configuration.)  Stress 
can  then  be  inferred  assuming  uniform  pressure  distribution 
over  the  known  sensors’  head  area. 

Sensors  are  mounted  on  the  surface  of  a  26  cm  diameter 
rigid  aluminum  wheel  (see  Figure  6).  Note  that  a  twin  wheel, 
without  the  array,  was  used  for  PIV  testing.  Five  sensors 
have  been  fabricated  and  integrated  in  a  linear  array 
spanning  one  half  of  the  wheel  width  (i.e.  from  one  edge  to 
the  center  of  the  wheel).  Sensors  were  first  calibrated  by 
applying  test  weights  of  100,  200,  and  500  grams  in  the 
normal  and  tangential  direction.  Measurement  linearity 
error,  across  all  the  sensors,  was  found  to  be  below  3%. 

The  sensor  array  is  extremely  sensitive  to  misalignment 
and  thus  an  uneven  contact  patch  profile  can  easily 
unbalance  the  output  reading.  To  ensure  accurate  alignment, 
sensors  alignment  was  verified  after  every  5  tests,  by  driving 
the  wheel  over  a  flat,  rigid,  aluminum  plate  covered  with  a 
thin  layer  of  polyurethane  foam  in  order  to  verify  that  the 
sensor  output  was  uniform.  Due  to  the  difficulty  in  precisely 
controlling  soil  preparation,  each  test  was  repeated  at  least 
15  times.  In  fact,  local  soil  density  variation,  inhomogeneity 
(due  to  non-uniform  distribution  of  larger  grains,  for 
example),  and  surface  unevenness  all  were  observed  to 
affect  measurement  output.  The  15  trials  highlighted  test 
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variability  and  were  analyzed  to  detect  outliers  and 
eventually  remove  tests  where  anomalies  were  detected. 

SOIL  PROPERTIES 

Characterization  of  the  soil  under  investigation  is  a 
necessary  step  for  any  terramechanics  investigation.  Detailed 
chemical  composition,  particle  size  distribution,  and 
shearing  properties  of  the  MMS  simulant  under  investigation 
can  be  found  in  [6].  However,  pressure-sinkage  properties 
(i.e.  Bekker’s  parameters)  for  the  soil  were  unknown,  and 
therefore  a  series  of  plate  penetration  tests  were  performed. 

Since  the  wheel  has  a  width  of  0.13  m  and  a  nominal 
contact  patch  length  of  0.05  m  (estimated  assuming  nominal 
conditions  of  Fz  =  100  N  and  low  slip)  three  rectangular 
plates  with  the  following  dimensions  were  selected:  0.13  m 
x  0.03  m,  0.13  m  x  0.05  m,  and  0.13  m  x  0.07  m. 

Each  plate  was  mounted  on  a  linear  actuator,  which  was 
anchored  to  the  testbed  and  then  pushed  perpendicularly  into 
the  soil  while  the  vertical  load  and  penetration  length  (i.e. 
sinkage)  were  measured  with  a  load  cell  and  a  draw  wire 
encoder,  respectively. 

For  each  plate,  tests  were  repeated  15  times.  Between  each 
test,  soil  was  manually  agitated  and  then  re-leveled.  Figure  7 
shows  an  example  of  the  data  collected.  Test-to-test 
variation  was  observed,  but  was  not  considered  unusual  due 
to  the  nondeterministic  nature  of  soil  testing. 

The  scope  of  the  tests  was  to  fit  experimental  data  to 
Bekker’s  pressure-sinkage  equation  [21]: 

P  =  (if  +  k<p)  zn  (3) 

where  p  is  pressure,  z  is  sinkage,  b  is  plate  width  (3,5,7  cm) 
and  {kc,  kphi,  n]  are  the  parameters  under  investigation. 
Adopting  the  fitting  methodology  presented  in  [2]  it  was 
noted  that  keqb  =  is  strongly  correlated  with  n  as 

shown  in  Figure  8.  This  correlation  necessarily  results  from 
the  tests  having  similar  amounts  of  deviation  from  an 
exponential  curve.  While  this  effect  is  solely  an  artifact  of 
experimental  estimation,  it  is  still  undesirable  because  it 
inhibits  keqb  from  being  estimated  independently. 

The  problem  is  mitigated  through  adoption  of  Reece’s 
equation  [22]  for  pressure-sinkage: 

P  =  keqr  (§)"  (4) 

Dimensional  analysis  of  Reece’s  equation  shows  that  keQr 
is  not  function  of  n  (as  it  was  in  Bekker  equation).  Although 
variability  is  still  substantial,  keQr  estimation  becomes  less 
dependent  of  n  as  can  be  seen  in  Figure  9. 
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Figure  7:  Penetration  tests  for  rectangular  plates  with  the  following 
dimensions  0.13  m  x  0.03  m,  0.13  m  x  0.05  m,  0.013  x  0.07  m. 


Figure  8:  Strong  correlation  between  soil  parameters  when  Bekker 
equation  is  used. 


Figure  9:  Correlation  between  soil  parameters  is  mitigated  when 
Reece's  equation  is  used. 
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Penetration  tests  variability,  even  under  laboratory 
controlled  conditions,  suggests  that  soil  parameters  should 
be  derived  from  statistical  distributions  rather  than 
deterministic  values.  A  stochastic  characterization  of  terrain 
properties  is  currently  being  investigated  by  the  authors 
while  the  results  presented  in  this  paper  are  still  derived  with 
the  method  established  by  Wong  [2]. 

Two  parameter  sets  are  reported  in  Table  1.  The  set  labeled 
‘357’  has  been  obtained  considering  the  full  dataset 
presented  in  Figure  7  while  the  set  labeled  “57”  has  been 
obtained  only  with  the  5  cm  and  7  cm  plates,  and  truncating 
the  data  at  50  kPa.  This  was  motivated  by  the  fact  that  the 
wheel  under  investigation  was  expected  to  have  contact 
patch  length  larger  than  5  cm  and  normal  stress  distribution 
below  50kPa.  The  two  datasets  show  how  slightly  modifying 
the  design  of  experiments,  can  drastically  change  soil 
parameter  calculation. 

Table  1 :  Bekker  soil  parameters  for  the  MMS  soil.  Two  sets  were 
extracted,  357  includes  all  the  data  while  for  57  only  two  plates 
were  used  (5  cm  and  7  cm)  and  data  was  truncated  at  50  kPa  mark. 

kc 

Set _ n _ [kN/mn+1]  [kN/mn+2] 

357  0.99  -55  4584 

57  1.4  846  6708 

RESULTS  AND  DISCUSSION 

Experiments  with  the  PIV  and  stress  sensor  experimental 
methodologies  were  conducted  separately.  For  PIV  tests,  a 
smooth  wheel,  coated  with  MMS  simulant  (to  ensure 
sufficient  interfacial  friction)  was  run  flush  against  a  glass 
wall.  For  stress  sensor  tests,  a  wheel  of  exactly  the  same 
diameter,  and  again  covered  with  MMS  simulant,  was  run  in 
the  middle  of  the  soil  bin.  Soil  was  loosened,  mixed,  and 
leveled  between  each  test,  in  an  attempt  to  achieve 
uniformly  loose,  homogenous  conditions. 

Both  type  of  tests  were  run  at  approximately  100N  of 
vertical  load  and  for  slip  levels  ranging  from  -70%  to  70% 
(for  PIV  tests,  slip  was  limited  to  ±30%).  For  PIV  tests  the 
wheel  velocity  was  fixed  at  17  deg/s  while  for  stress  sensor 
tests  angular  velocity  was  reduced  to  8.5  deg/s  to  improve 
measurement  quality.  (The  horizontal  carriage  velocity  was 
modified  to  achieve  the  desired  slip  level.)  For  both  types  of 
tests,  it  was  first  ascertained  that  velocity  did  not  have  an 
influence  on  wheel  performance.  The  operational  conditions 
described  above  were  chosen  because  they  are  close  to  those 
of  the  Mars  Exploration  Rover,  a  successful  lightweight 
robotic  vehicle. 

A  substantial  amount  of  data  was  collected  and  cannot  be 
comprehensively  described  in  this  paper.  Instead,  a  small 
number  of  initial  results  are  presented. 


Figure  10:  A  snapshot  of  a  30%  slip  test.  Nominal  vertical  load 
was  100N  and  wheel  angular  velocity  of  17  deg/s.  From  top-left- 
clockwise:  velocity  vectors,  u- velocity,  v- velocity,  and  velocity 
magnitude. 


PIV  Analysis 

Analysis  of  PIV  data  was  performed  to  qualitatively 
analyze  soil  motion  (a  quantitative  analysis  would  have 
required  to  investigate  the  complex  mapping  between  stress 
and  displacement,  this  goes  beyond  the  scope  of  this 
preliminary  study).  Figure  10  presents  a  snapshot  of  a  30% 
slip  test,  and  displays  the  following  information  from  top- 
left-clockwise:  velocity  vectors,  u- velocity,  v- velocity,  and 
velocity  magnitude.  Analysis  of  such  images  can  provide 
insights  into  the  spatial  distribution  of  soil  velocity  under 
running  gear,  and  can  vary  dramatically  for  such  cases  as 
slip,  skid,  free-rolling  wheels,  braked  wheels,  etc. 

Decomposition  of  this  flow  field  can  yield  useful  insight 
into  soil  shearing  (which  occurs  primarily  in  the  horizontal 
direction,  see  upper  right  image)  and  soil  compaction 
phenomena  (which  occurs  primarily  in  the  vertical  direction, 
see  lower  right  image).  Here,  a  blue  region  corresponds  to 
no  motion  while  red  indicates  a  maximum  velocity.  Analysis 
of  these  images  shows  that  soil  flow  remains  attached  to  the 
wheel  rim.  Moreover,  for  low  vertical  load  (such  as  the  one 
utilized  during  experiments)  it  was  observed  that  two 
separate  slip  failure  lines  did  not  evolve,  as  predicted  by 
classical  theory  [23,  24].  This  finding  is  interesting  because 
according  to  [23],  the  maximum  stress  occurs  where  the  soil 
flow  separates.  The  absence  of  flow  separation,  however, 
does  not  prevent  stress  to  reach  a  maximum  (see  Figure  11). 

For  slip  levels  below  ±10%,  the  soil  was  not  observed  to 
develop  a  significant  shearing  plane.  Another  phenomenon 
that  was  clearly  highlighted  by  PIV  analysis  is  the  periodic 
nature  of  soil  failure.  For  slip  level  above  10-15%,  soil  often 
exhibits  a  periodic  loading  cycle  of  alternating  compaction 
and  shearing,  which  results  in  discontinuous  failure  of  the 
soil  mass.  This  has  two  direct  consequences:  oscillations  in 
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drawbar  pull  readings  and  creation  of  ripples  behind  the 
wheel.  Note  that  while  these  effects  have  been  noted 
previously,  they  have  been  typically  assigned  to  the  effect  of 
grousers.  However,  these  effects  are  present  even  for  smooth 
wheels,  without  grousers. 

PIV  data  can  be  useful  for  investigation  of  constitutive 
models  for  granular  materials,  and  for  development  of 
reduced  order  models  based  on  soil  displacement 
predictions.  An  important  consideration  to  bear  in  mind 
when  examining  flow  fields  like  the  one  presented  in  Figure 
10  is  that  the  relationship  between  stress  and  displacement  is 
typically  complex,  and  one  must  avoid  the  temptation  to 


directly  (i.e.,  proportionally)  correlate  velocity  magnitudes 
with  stress  magnitudes. 

For  this  reason,  direct  stress  measurement  of  shear  and 
normal  forces,  and  inferences  of  associated  stresses,  at  the 
wheel-terrain  interface  yields  valuable  information  about  the 
traction  generation  process. 

Interface  Force  Sensor  Analysis 

Classical  terramechanics  methods  rely  on  the  estimation  of 
the  stress  distribution  under  the  wheel.  The  ability  to  directly 
measure  such  quantities  allows  for  a  one-to-one  comparison 
of  model  prediction  and  experimental  reality. 


Figure  11:  Normal  and  tangential  stress  at  the  wheel-soil  interface  calculated  from  force  sensors.  These  were  obtained  for  nominal  vertical  load 
of  100  N  and  wheel  angular  velocity  of  8.5  deg/s.  The  four  panels  present  data  for  -70%,  +70%,  -10%,  and  +10%  slip  (clockwise  from  upper 
left).  Sensors  are  labeled  according  to  the  scheme  presented  in  Figure  6.  “I”  corresponds  to  the  sensor  located  at  the  center  and  “V”  to  the  sensor 
located  at  the  edge  of  the  wheel.  Central  angle  defines  the  angular  position  along  wheel  circumference  [26,27]. 
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Analysis  of  stress  distribution  across  a  (symmetric)  half- 
wheel  width  shows  that  boundary  effects  become  more 
pronounced  as  slip  increases  (see  Figure  11).  In  particular, 
stress  at  the  wheel  edge  was  observed  to  be  relatively  high 
for  positive  slip  and  relatively  low  for  negative  slip.  It  is 
hypothesized  that  this  effect  is  caused  by  soil  transport 
phenomena:  for  positive  slip,  soil  in  the  center  of  the  wheel 
is  transported  behind  the  wheel  at  higher  rate  than  the  soil  at 
the  wheel  edges,  which  causes  the  wheel  edges  to  bear 
proportionally  more  of  the  total  normal  wheel  load.  On  the 
contrary,  for  negative  slip,  soil  accumulated  in  front  of  the 
wheel  creates  a  thicker  layer  under  the  wheel  median  axis, 
causing  higher  stress  in  the  center. 

For  higher  loading  conditions,  Onafeko  and  Reece  [25] 
noted  that  normal  stress  decreases  with  increasing  positive 
slip  since  an  increasingly  larger  portion  of  vertical  load  is 
supported  by  shear  stress  (which  contributes  more  to  vertical 
load  equilibrium  because  of  increased  sinkage).  This  was 
confirmed  experimentally  with  the  stress  sensors. 

Another  interesting  aspect  of  wheel  stress  distributions  is 
the  inversion  of  shear  stress  for  negative  slip  conditions. 
This  phenomenon  was  noted  also  by  [25]  and  it  is  consistent 
with  wheel-soil  interaction  kinematics:  for  negative  slip,  the 
wheel  travels  forward  but  simultaneously  skids  over  the  soil, 
generating  a  shear  sign  transition.  Interestingly,  PIV  imagery 
does  not  show  any  soil  separation  or  flow  inversion  where 
the  shear  stress  changes  sign. 

In  Figure  12,  a  direct  comparison  between  the  measured 
stress  and  stress  predicted  by  the  model  originally  proposed 


by  Wong  [26,  27]  and  Janosi  and  Hanamoto  [28]  is 
presented,  using  the  experimentally  determined  soil 
parameters  (two  parameter  sets,  presented  in  Table  1,  are 
compared).  The  normal  stress  distribution  is  underestimated 
and  the  error  seems  largely  related  with  the  location  of 
maximum  stress.  Tuning  of  semi-empirical  model 
parameters  could  allow  better  agreement. 

The  predicted  shear  stress,  however,  was  found  to  be 
overestimated.  Note  that  the  shear  modulus  adopted  to 
produce  results  in  Figure  12  was  calculated  according  to 
[29].  For  larger  (but  arguably  inaccurate)  values  of  shear 
modulus,  it  may  be  possible  to  obtain  better  agreement 
between  prediction  and  experimental  data;  however  this 
raises  a  fundamental  question  about  the  validity  of  the 
assumptions  behind  the  model.  In  fact,  the  model  assumes 
that  the  soil  is  sheared  for  a  distance  corresponding  to  the 
amount  of  relative  motion  between  the  wheel  and  the  soil. 
This  assumption,  as  shown  by  PIV  analysis,  is  likely 
erroneous,  since  the  soil  at  the  wheel-terrain  interface  stays 
attached  to  wheel  rim,  while  failure  physically  occurs  (in 
regular,  periodic  failure  patterns)  some  distance  away  from 
the  interface.  Although  [n,  kc,  k^}  and  {n, /tc,k0}57  are 
significantly  different  (see  Table  1),  model  predictions  using 
these  two  sets  are  relatively  close.  This  warrants  further 
efforts  in  characterizing  terrain  variability  and  its  influence 
on  stress  measurements  variability. 


Figure  12:  Stress  distribution  for  10%  (left)  and  30%  (right)  slip  compared  with  analytical  model  from  Wong  and  Reece  [26,  27].  Two  soil 
parameter  sets,  presented  in  Table  1,  were  tested.  The  difference  between  the  two  parameter  sets,  although  significant,  it  is  not  dramatic.  Normal 
stress  is  slightly  underestimated  while  tangential  stress  is  significantly  estimated.  Tangential  stress,  however,  is  primarily  based  on  soil  shear 
properties  which  were  obtained  in  [29]. 
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CONCLUSIONS 

Novel  experimental  methods  aimed  at  understanding  the 
fundamental  phenomena  governing  the  motion  of 
lightweight  vehicles  on  dry,  granular  soils  were  presented. 

Aside  from  standard  wheel  experiments  (i.e., 
measurements  of  drawbar  force,  applied  torque,  and  sinkage 
during  controlled  slip  runs)  two  additional  experimental 
methodologies  were  introduced.  The  first  relies  on  high¬ 
speed  imaging  of  the  wheel- soil  interface  and  the  use  of 
particle  image  velocimetry  (PIV)  to  measure  micro-scale 
terrain  kinematics.  The  second  experimental  methodology 
consists  of  a  custom  force  sensor  array  located  at  the  wheel- 
terrain  interface.  The  sensors  allowed  explicit  measurement 
of  normal  and  shear  forces  (and,  therefore,  estimation  of 
normal  and  shear  stresses)  at  numerous  discrete  points  along 
the  wheel-soil  interface. 

Analysis  of  PIV  data  has  shown  that  soil  failure,  at  certain 
slip  levels,  is  qualitatively  different  under  cases  of  low 
vertical  load  (which  is  typical  for  lightweight  robots) 
compared  to  cases  of  high  vertical  load  (typical  for  large 
ground  vehicles).  Also,  soil  flow  patterns  have  been 
observed  to  exhibit  periodic  failure  phenomena,  giving  rise 
to  interesting  features  such  as  surface  ripple  formation.  Soil 
flow  was  observed  to  be  always  attached  to  the  wheel  rim 
and  only  one  shear  failure  surface  was  observed.  Soil  usually 
exhibits  compression  in  front  of  the  wheel  and  then  shears 
beneath  it. 

Stress  measurements  showed  that,  although  only  one  shear 
failure  surface  is  present,  tangential  stress  goes  through  sign 
inversion  for  negative  slip.  Stress  distribution,  along  the 
wheel  width,  is  approximately  uniform  for  low  slip  while 
edge  effects  become  increasingly  significant  for  higher  slip 
levels.  Although  some  observations  regarding  soil  shear 
failure  were  not  confirmed  by  PIV,  classical  methods 
(partially  based  on  those  observations)  were  able  to  capture 
main  trends  for  a  range  of  slip  conditions.  These  results 
provide  deeper  understanding  of  the  mechanics  of  traction 
generation  and  are  expected  to  open  new  frontiers  for  more 
accurate,  and  predictive,  lightweight  vehicle  mobility 
models. 

Further  investigation  of  small  robot-terrain  interaction 
mechanics  will  focus  on  extending  these  experiments  to  a 
wider  range  of  vertical  loads.  This  will  provide  a  basis  for 
validation  of  constitutive  laws  and  the  improvement  of 
reduced-order  models.  Future  work  will  also  focus  on 
stochastic  characterization  of  terrain  response  and  how 
underlying  soil  variability  affects  interfacial  stresses 
modeling.  In  fact,  even  under  laboratory  controlled 
conditions,  penetration  plate  tests  have  highlighted 
significant  soil  variability,  warranting  for  statistical 
interpretation  of  experimental  data. 
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